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Future large scale structure observations are expected to be sensitive to small neutrino masses, of 
the order of 0.05 eV or more. However, forecasts are based on the assumption that by the time 
at which these datasets will be available, the non-linear spectrum in presence of neutrino mass 
will be predicted with an accuracy at least equal to the neutrino mass effect itself, i.e. about 3%. 
Motivated by these considerations, we present the computation of the non-linear power spectrum of 
ACDM models in the presence of massive neutrinos using the Renormalization Group (RG) time- 
flow approach, which amounts to a resummation of perturbative corrections to the matter power 
spectrum to all orders. We compare our results with those obtained with other methods, i.e. linear 
theory, one-loop perturbation theory and N-body simulations and show that the time-RG method 
improves the one-loop method in fitting the N-body data, especially in determining the suppression 
of the matter power spectrum when neutrino are massive with respect to the linear power spectrum. 

PACS numbers: 98.80.Cq 



I. INTRODUCTION 

The detection of atmospheric neutrino oscillations proves that at least one neutrino mass eigenstate is as heavy as 
approximately (2.4 ± 0.1 x 10~^eV^)^/^. In the inverted hierarchy scenario, two eigenstates would have such a large 
mass, while in the degenerate scenario, the three masses would be of this order of magnitude or larger (for a recent 
analysis of the data, see e.g. [1]). In all cases, the total mass M^, = J2 "^i^ should be greater than approximately 0.05 
eV. This means that neutrino masses are expected to impact significantly the cosmological evolution during matter and 
dark energy domination, first by enhancing the background neutrino density (today, ^li, = {M^/93AAeV)h~^ ^ 10""^, 
instead of fi^ = 1.712 x 10~^h~^ in the massless limit) [2J, and more importantly, by reducing the growth rate of 
matter perturbations on comoving wavenumbers larger than /cm. ~ 2 x lO^'^ft./Mpc. In the linear theory, this effect is 
maximal for k <^ 1/i/Mpc and suppresses the current matter power spectrum (PS) by a factor (1 — 8^1^ i-e. at 
least by 3% [3, 4,. 

At the moment, the most robust cosmological bounds on the total neutrino mass M^, — '^m^ rely either on its 
background effect (Komatsu et al. obtain < 0.67eV [95%C.L.] with WMAP-5, BAO and SNIa data only), or on 
limits on the suppression of the matter PS in the linear regime (Tegmark et al. |S] obtain M,^ < 0.9eV [95%C.L.] using 
the combination of WMAP-3 with the PS of SDSS-LRG in the range 0.01 < k < 0.2/i/Mpc, after marginalization over 
an unknown nuisance parameter describing non-linear effects in the range 0.1 < k < 0.2/i/Mpc). In order to improve 
significantly these limits, it would be crucial to model the effect of neutrino masses in the non- linear regime, in order to 
include smaller scales in the analysis and increase the lever arm. For galaxy redshift surveys, including wavenumbers 
above k 0.1/i/Mpc remains a delicate issue, because in addition to non-linear corrections to the theoretical matter 
PS, one should be able to predict the scale-dependence of the light-to-mass bias, which starts to play a role on those 
scales. However, the prospects for neutrino mass determination are particularly encouraging thanks to other types of 
experiments. For instance, future weak lensing surveys (combined eventually with CMB lensing measurements) might 
reach a sensitivity of Mi, 0.04eV (2-cr) according to [7] (see also [HI IHl IIHI HI])- These measurements - which are 
not affected by any light-to-mass bias - can probe the matter PS up to z ~ 3, which means, first, that they are likely 
to detect the non-trivial redshift dependence of the small-scale matter PS in the presence of massive neutrinos; and 
second, that they can measure cosmological fluctuations deeper in the linear regime than observations limited to small 
redshifts. Nevertheless, the forecasts of [TIIHIIS] do include non-linear scales, corresponding to multipoles £ < 1000 in 
the case of [7] (for comparision, non-linear effects are important a.t £ > AO for sources with Zs ^ 0.2, and £ > 200 for 

~ 3). Hence, they are based on the assumption that by the time at which these datasets will be available (of the 
order of one decade) , the non- linear spectrum in the presence of neutrino mass will be predicted with an accuracy at 
least equal to the neutrino mass effect itself, i.e. about 3%. The same comment applies to other promising forecasts 
related to future cluster surveys [12], 21cm surveys [131 IHIj high-redshift galaxy surveys [IF, CMB-weak lensing 
cross-correlation ^TBJITT] or Lyman-a forest data [TB]. In the latter case, one should model the effect of neutrino mass 
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not only on the total matter power spectrum, but also on the thermodynamical evolution of the intergalactic medium: 
hence, massive neutrinos should be implemented in hydrodynamical simulations of structure formation before deriving 
robust bounds on neutrino masses. All previous neutrino mass bounds from Lyman-a forest data neglected this issue. 

Current weak lensing data probe the matter PS on much smaller scales than the previously discussed experiments 
(typically of the order of A; ~ (10 — 40) ft./Mpc). In this case, neutrino mass bounds have been obtained using the 
HALOFIT approach in order to compute non-linear corrections to the CDM component of the total power spectrum 
[iniHn]- The HALOFIT algorithm was derived and calibrated for the ACDM model. Its use in presence of massive 
neutrinos might be correct, but this would need further justification (see |21j). 

Finally, the computation of the non-linear PS in presence of massive neutrinos would be important for CMB 
physics. A crucial issue in CMB observations is to deal with foreground contamination. In particular, the thermal 
SZ effect is difficult to remove, and a large fraction of it will contribute to the final measured temperature spectrum. 
This SZ effect depends on the non-linear PS; hence, it can be computed and fitted at the same time as primary 
anisotropies. However, any mistake in the computation of this contribution would result in a significant bias on the 
cosmological parameters (for instance, on the scalar spectral tilt) |22]. Even with the smallest allowed neutrino mass, 
the SZ spectrum should be sensitive to neutrino free-streaming on small scales. So, if the CMB spectrum of a precise 
experiment like Planck was analysed using SZ spectrum computations in which neutrino mass effects are neglected, 
the bias induced on cosmological parameters might be significant. Note that instead of fitting this component, one 
can marginalize over nuisance parameters accounting for its amplitude and slope, or just restrict the multipole range 
used in the analysis, but this is at the expense of reducing the parameter sensitivity |22j . Hence, investigating the 
effect of neutrino masses on non-linear cosmological perturbations is also important in order to push the precision of 
CMB observations on all comsmological parameters. 

For all these reasons, a crucial step for precision cosmology is to incorporate free-streaming massive neutrinos in 
numerical or analytical estimates of the non-linear matter PS on small scales. On the side of numerical simulations, 
the first steps were performed in [331 [211 US . To our knowledge, other attempts to simulate mixed cold plus hot 
models were treating neutrinos as extra cold particles with a different initial PS. This approach neglects the effect of 
neutrino thermal velocities, and fails to reproduce even the linear evolution. Ref. [25] shows explicitely that neglecting 
thermal velocities results in an error which is of the same size as the very effect of a small neutrino mass: this makes 
current bounds on nii, from Lyman-a data suspicious. However, implementing free-streaming neutrino "particles" in 
simulation is very demanding, because of the large number of particles required to sample the full neutrino phase 
space, and of the tiny time steps needed to follow the particle trajectories with sufficient accuracy, because of their 
large velocities. Some significant progress can be expected in this field in the next years. 

In parallel, it is important to derive analytic approximations of the non-linear evolution in order to cross-check the 
results from simulations, to reach eventually smaller scales and larger redhifts, and hopefully to obtain computationally 
faster ways of predicting non-linear spectra. Various analytical schemes have been discussed in the past years for 
ACDM models, inspired by quantum field theory; these include one-loop calculations [371 UHl HHj) resummation 
techniques [30l [311 132], renormalization group aproaches [33l O |35l |36l |37l |38l |39|, or the Time Renormalization 
Group (TRG) flow calculation of [40^. 

The first attempt to incoporate massive neutrinos in these calculations was performed by Saito et al. ^41j . who 
provide a first-order approximation of the one-loop computation. The full one-loop results were obtained by Wong 
[42j . The goal of this paper is to generalize the method of [40' in the presence of neutrino masses, and to compare with 
previous results from [25, 41, 42J. In this new approach the matter power spectrum, the bispectrum (i.e. the connected 
three point-correlation function in Fourier space) and higher order correlations, are obtained - at any redshift and for 
any momentum scale - by integrating a system of differential equations. Truncating at the level of the trispectrum 
(i.e. the connected four point-correlation function in Fourier space), the solution of the equations corresponds to 
the summation of an infinite class of perturbative corrections. In the diagrammatic language, this is equivalent 
to resumming at all orders in perturbation theory the corrections to the matter power spectrum. Furthermore, 
comparison with N-body simulations for zero neutrino masses shows that the method is valid up to momentum 
wavenumbers of the order of 0.3 /i/Mpc at z = 1, and able to describe the non-linearities due to mode-to-mode 
coupling in the baryon acoustic oscillations region. Compared to other resummation frameworks, this scheme is 
particularly suited to cosmologies other than ACDM, such as those based on modifications of gravity and those 
containing massive neutrinos. 

The paper is organized as follows: in Section 1 we summarize our method based on Ref. |40j generalized to include 
massive neutrinos. In Section 3 we present our numerical findings and compare them with those obtained with other 
methods, namely linear theory, one-loop perturbation theory [311133] and the N-body simulations of Refs. [26J. 
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II. METHOD 



Our goal is to set up a method to compute non-linear corrections to the matter PS in a mixed dark matter (MDM) 
model. The matter density fluctuation is given by 

<5m = {Spc + 6pb + ^Pv)lp,n = (1 - !v)^cb + ! , (1) 

where the subscripts 'c', '6', V, and 'c6', stand for cold dark matter (CDM), baryons, neutrinos, and CDM plus 
baryons, respectively, and f^, = Cl^/flm- The total matter PS is then given by 

Pm{k; z) = {l- Uf P,i,{k; z) + 2 (1 - U)U PcbAk; z) + fl Pu{k; z). (2) 

In ref. |41j the expression above was approximated by using the linear order result for the cb — v and the v — v PS's 
(the second and third terms in Eq. ([2| , respectively) , while including non-linear corrections to the ch — cb one in the 
so-called one-loop approximation, which includes density and velocity perturbations up to the third order |27j . In 
the present paper we will go beyond the one-loop approximation for Pcb , using the TRG approach of Ref. |40j , while 
keeping the linear order results for the other two contributions to Eq. 

The TRG is a method to sum perturbative corrections to all orders. The starting point are the continuity and 
Euler equations satisfied by the density contrast and peculiar velocity of a non-relativistic fluid {i.e. CDM or baryons 
in the present case), 

?^+V.[(H-<5,fc)v] = 0, 

O T 

-^ + 'Wv +(v V)v= -V0, (3) 
a T 

where r is the conformal time and we have assumed db = Sc = 6cb. The gravitational potential </> is determined by the 
total mass fiuctuation, via the Poisson equation 

v^^^^n^n^s^, (4) 

where Ti. — dlog a/ dr . Dropping the ^cfcv term in the continuity equation and the (v • V)v one in the Euler equation, 
the solutions to the system above reproduce linear perturbation theory on the subhorizon scales we are interested in. 

In order to close the system ([3|, we go to Fourier space, use Eq. ([T|, and approximate the RHS of the Poisson 
equation as follows 



^ n^{T) S^{k, t) ~^H^ nll^ik, t) <5,,(k, r) , (5) 



where 



(1 - .U)Sg{K r) 



nli^ik, r) ^ f^„(r)(l -/.)(!+ ""7^71!.' ^ ) • (6) 



(5^^^ indicate the density perturbations evolved according to linear theory. Notice that, due to the different space 

dependence of S^; and 5"^ in the massive neutrino case, - unlike Q.m - is space-dependent. Another approximation 
is to recursively use the full TRG-evolved S^b fluctuation in Eq. Q. However, we numerically checked that performing 
one iteration (i.e., using the output non-linear 5cb as an input in Eq. (l6|) instead of 5^b) does not affect the results by 
more than 0.1%. This way of including massive neutrinos is conceptually similar to the grid method used in N-body 
simulations of Ref. [55]. where, instead of simulating massive neutrinos as particles with individual velocities, they 
are embedded as a local neutrino density on a grid which is evolved in time using linear theory. 

Next, we write Eqs. (|3] |4]) in a compact form [Snilin]- First, we introduce the doublet Lpa (a = 1,2), given by 

¥'i(k,??) \ = ^-n ( Scb{k,T]) 



— ^ { -Mir JiW-w ' (''') 



where = ik • v, and the time variable has been replaced by the logarithm of the scale factor. 



r? = log — , (8) 
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ain being the scale factor at a conveniently remote epoch, such that all the relevant scales are well inside the linear 
regime. 



Then, we get 



where 




and the only non- vanishing elements of the vertex function 7£if,c(k, p, q) are 

7ii2(k, p, q = - k + p + q ^ 

2 q"^ 

{p + qfp- q 



7222(k, p, q) ^ Soik + p + q) f ^ , (11) 



and 7121 (k, p, q) = 7112 (k, q, p). 

In Eq. (|9|, repeated indices are summed over, and integration over momenta q and p is understood. Notice that 
all the information on the neutrino mass, both at the background and at the linear perturbation level, is contained 
in O21, the other entries of and the vertices being universal. 

At this point we should stress a crucial difference between the present approach and the ones of refs. [30l [3TJ [Ml [35] . 
Namely, in those papers, the equations are derived in the Einstein-deSitter (EdS) case = 1) said then extended to 
different cosmologies, such as ACDM, by reinterpreting i] as the logarithm of the linear growth factor, while keeping 
the equations unchanged. In particular, the - matrix is approximated by that of the EdS model. It was discussed 
in [30' and numerically checked in |40j that in cosmologies with a constant equation of state for the dark energy, this 
procedure gives at most a (9(1%) error at z = for A; ^ 0.3ft./Mpc, rapidly decreasing at higher redshift and larger 
scales. The physical reason for this accuracy lies in the fact that this procedure takes fully into account the growing 
mode of the non-EdS cosmology, while it mistreats the decreasing one. The latter start to play a role only at high k's 
or for low redshifts, where non-linearities become important. 

However, when massive neutrinos contribute to the dark matter, the linear growth factor is fc-dependent, and the 
redefinition of rj is ill-defined. This singles out the present approach, the TRG, as particularly suited to this case. 
Indeed, in this approach, we keep the same definition, Eq. ((sl) for any cosmological model, and take fully into account 



the fc-dependence of the linear growth function via the 51^^ (k, 77) term in Eq. (10). 

The evolution in rj of the correlation functions of the ipa fields can be derived by iterating the application of Eq. (|9| . 
The result is an infinite tower of coupled integro-differential equations. Following Ref. [40], we will truncate the system 
by neglecting the trispectrum. The evolution equations for the PS, (iy9a(k, '7)vfc(qj v)) = ^d(M + q)^a6(k,?7), and for 
the bispectrum, {ipaO^, v)Vb{(l, v)'Pc{p, v)) = '^^(k + q + p)Babc{k, q, p; v) are then given by 

+ J [7acd(k, -q, q - k) Bbcdi^, -q, q - k; 77) 
+ Bacdi'k, -q, q - k; 77) 7hcd(k, --q, q - k)] , 



dj^ BabciK -q, q - k; 77) = -6ad(k , ifjBdbdK -q, q - k; ry) - 66d(-q , v)BadciK -q, q - k; 77) 

-0c<i(q - k , r])Babdik, -q, q - k; 77) 
+2e'' [7ade(k, -q, q - k)P<jb(q , v)Peci'k - q , ??) 
+76de(-q, q-k, k) Pdc (k - q,7y) Pea (k, 77) 

-H7cde(q-k, k, -q)P<j„(k,7?)Pef,(q,7?)] . (12) 

The system is solved by giving initial conditions for the PS and the bispectrum at 77 = 0, and integrating forward in 
?/. In the ^ limit, the initial conditions can be obtained from linear perturbation theory, 

Pii(k, 0) = P/l(k, 0), 
Pi2(k, 0) = /L(k,0)P/l(k, 0), 
P22(k, 0) = /|(k,0)P/l(k, 0), 

P„6c(k, -q, q-k; 0) =0, (13) 
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FIG. 1: Matter power spectra at a = 0.3 {z = 2.33) normalized to the linear one for Mi, = 0. The red, solid lines are the result 
of the TRG method described in the text, while the green, dash-dotted ones are the linear approximation. The curves from 
top to bottom correspond to = 0, 0.3, and 0.6 eV, respectively. 



where the multiplication of the initial P12 and P22 by appropriate powers of the linear growth function, /L(k, 77) = 
dlog(5^j,(k, ri)/dri, selects the linear growing mode. In our computations, we will use = 1/41, and take the density- 
density power spectrum from the Boltzman code CAMB J43]. More precisely, we extract Pcb, Pcb,v, Pv and ri^^^ 
at any value of (fc, z) from the quantities 5cb and Sv evolved by CAMB in the synchronous gauge. 

Different approximations of the system above give the linear and the one- loop approximations, respectively. Indeed, 
the linear result for the PS is obtained simply by setting Babe = 0, and integrating the first of Eqs. ( 12 1. On the other 
hand, the one-loop result is obtained by inserting the linear PS in the RHS of the second of Eqs. 
resulting Babe in the first one. 

Integrating the system above as it is, goes beyond the one-loop approximation, and corresponds to summing infinite 
perturbative contributions. As discussed in detail in Ref. |40| . at each step in 77 the PS and the bispectrum appearing 
at the RHS's include non-linear corrections maturated from 77 = to ry. 

In the next section we will present the results of such integrations and compare them with the linear, the one-loop 
approximations and the N-body results. 



12 1 and using the 



III. RESULTS AND DISCUSSION 



The impact of non- linear corrections on the linear PS is opposite with respect to that of a massive neutrino component. 
Indeed, while the latter damps the PS at scales larger than ~ 2 x 10~^/i/Mpc, non-linearities enhance it for 
k ^ A:nl, where /cnl — 0.05/i/Mpc at z = and scales roughly as (1 -I- z)^/^. The combination of both effects is 
shown in Fig. [T] where we plot the total matter PS normalized to the linear one for M^, = 0, at a = 0.3, i.e. z — 2.33 
(a is the cosmological scale factor, normalized to a = 1 today). We have assumed a flat cosmological model with 
Hubble parameter h — 0.70, and density parameters — 0.7 and fi^ — 0.3. fi^ is given by the sum of ^b = 0.05 
and + f^iy. The primordial power spectrum was assumed of the scale-invariant Harrison-Zel'dovich form (n = 1), 
normalized as to give cts = 0.878 for Q.^ = 0. 

The green, dash-dotted lines are the linear approximation result, while the red ones were obtained with the TRG 
method. The curves from top to bottom correspond to — 0, 0.3, and 0.6 eV, respectively. As we see, even at 
an intermediate redshift the effect of non-linearities overcomes that of massive neutrinos, turning a fc-independent 
suppression in a /c-dependent enhancement, which should be understood in order to derive reliable bounds on neutrino 
masses. 

In the following figures we disentangle the non-linear from the massive neutrinos effects. Non-linear effects are 
shown in Fig. |2] where the total matter power spectra is divided by the linear one computed for the same value of the 
total neutrino mass. We compare results for M^, = 0, 0.3 and 0.6 eV obtained in one-loop perturbation theory (black 
dotted) , and by the TRG method (red solid) . Blue diamonds are the results of N-body simulations of Ref. [25] , when 
neutrinos are treated like particles. The results are given for a = 0.3 corresponding to redshift z = 2.33 and a = 0.5 
corresponding to z = 1. 
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FIG. 2: Power spectra divided by the linear one for the same value of Mv as computed in one-loop perturbation theory (black, 
dotted), and by the TRG method (red, solid). Blue diamonds are the results of N-body simulation of Ref. (26j, when neutrinos 
are treated like particles. 



At high redshift, the TRG method works surprisingly well in the whole range of k displayed in the figure, and is 
able to reproduce the N-body simulations better than the one-loop result, especially for large neutrino masses. For 
smaller redshifts, this tendency remains, but as expected, the agreement between the N-body simulation and the TRG 
method is valid only up to fc = 0.3/i/Mpc, while in the same range the one-loop method overestimates the non-linear 
to linear power spectra ratio. At both redshifts, the agreement with N-body simulations increases with the neutrino 
mass. This can be understood by recalling that we keep fie + fixed, and non-linear effects decrease when the CDM 
component is lowered. 

In Fig.|3]we show the effect of the neutrino mass as computed in linear approximation (green, dash-dotted line), one- 
loop perturbation theory (black, dotted), and by the TRG method (red, solid). Blue diamonds are again the results 
of N-body simulation of Ref. [26j.The upper (lower) group of lines corresponds to A/^ — 0.3 (0.6) eV, respectively. 
From this figure, we can infer two things. On one hand, the relative change between massive and massless neutrino 
models is approximately the same for the TRG and the one-loop method. This is not surprising, because the difficulty 
of the one-loop method to reproduce the N-body data disappears in the plotted ratios. On the other hand, we see 
that both the one-loop and RG method tend to over-estimate the impact of neutrino masses for small redshift. It 
should also be noticed that the simulations of Ref. [26], where the neutrinos are modeled by the grid approach, seem 
to be in closer agreement with us than those in which neutrinos are represented as particles (the blue diamonds in 
our plots), as could be expected from the conceptual similarity. 

In this paper we have presented the derivation of the large scale matter power spectrum in the presence of massive 
neutrinos using the TRG method which resums the perturbative corrections to the matter power spectrum to all 
orders and therefore extends the one-loop computation. Evaluating numerically the TRG matter power spectrum, 
we have seen that nonlinear corrections enhance the suppression of the power spectrum at small scales due to the 
neutrino free-streeming. This enhanced suppression has been observed in N-body simulations |25} I26j and already 
appears at one-loop [42 , but the TRG is in better agreement with N-body data, as it is manifests from Fig. 2. 

Our results indicate that the TRG approach is particularly suited to describe intermediate scales, where linear 
theory and one-loop perturbation theory are not reliable, and, on the other hand, N-body simulations require very 
large volumes and are therefore extremely time consuming. This range of scales will be accurately measured by 
future galaxy surveys with the aim of pinning down the dark energy equation of state by measuring baryonic acoustic 
oscillations (BAO), from k ~ 0.05/i/Mpc to 0.3/i/Mpc. In this context, it has been shown [3T1I33] that non-linearities 
affect non-trivially the BAO feature by shifting the position of the peaks towards smaller scales with respect to those 
computed in linear theory. The same effect is witnessed by the wiggles in Figs.[2j which would not be there if the non- 
linear corrections were /c-independent in that range of scales. Neutrino masses, modifiying the non-linear corrections, 
also have an impact on the BAOs, and this is seen in our Figs. |3] 
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FIG. 3: Suppression on the total PS induced by non-zero neutrino masses at difTerent epochs, as computed in linear approx- 
imation (green, dash-dotted line), one-loop perturbation theory (black, dotted), and by the TRG method (red, solid). Blue 
diamonds are the results of N-body simulation of Ref. [5^. The upper (lower) group of lines corresponds to Mi, — 0.3 (0.6) eV, 
respectively. 



The approach can be improved along two possible lines. On one side neutrino perturbations have been treated at 
the linear level; a nonlinear generalization would be of course welcome. On the other side, the TRG technique might 
be improved by including vertex corrections or, in other words, solving the equations for the correlators including the 
connected four-point correlation function. 
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